Multivariate analysis (bam)
Logistic regression models
Birthweight
#baseline model
m_BW1 <- lm(BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex, data=bevn_eco_in7)
summary(m_BW1)##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex,
## data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2701.0 -264.2 -13.3 249.7 5614.7
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.717e+03 8.853e+00 -419.88 <2e-16 ***
## GA_weeks 1.790e+02 2.237e-01 800.35 <2e-16 ***
## birth_Y_M_num -1.060e-01 7.876e-03 -13.46 <2e-16 ***
## parity_cat(1,2] 1.349e+02 8.198e-01 164.62 <2e-16 ***
## parity_cat(2,3] 1.787e+02 1.253e+00 142.59 <2e-16 ***
## parity_cat(3,20] 1.991e+02 2.189e+00 90.97 <2e-16 ***
## sex2 -1.470e+02 7.515e-01 -195.55 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 394.4 on 1103345 degrees of freedom
## Multiple R-squared: 0.3895, Adjusted R-squared: 0.3895
## F-statistic: 1.173e+05 on 6 and 1103345 DF, p-value: < 2.2e-16 #with mat age
m_BW1_matage_cont = update(m_BW1, . ~ . + mat_age)
summary(m_BW1_matage_cont)##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex +
## mat_age, data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2701.1 -264.2 -13.3 249.8 5613.8
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.712e+03 9.284e+00 -399.810 <2e-16 ***
## GA_weeks 1.790e+02 2.241e-01 798.992 <2e-16 ***
## birth_Y_M_num -1.048e-01 7.896e-03 -13.276 <2e-16 ***
## parity_cat(1,2] 1.352e+02 8.324e-01 162.471 <2e-16 ***
## parity_cat(2,3] 1.792e+02 1.276e+00 140.428 <2e-16 ***
## parity_cat(3,20] 1.999e+02 2.219e+00 90.078 <2e-16 ***
## sex2 -1.469e+02 7.515e-01 -195.544 <2e-16 ***
## mat_age -1.582e-01 7.788e-02 -2.031 0.0422 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 394.4 on 1103344 degrees of freedom
## Multiple R-squared: 0.3895, Adjusted R-squared: 0.3895
## F-statistic: 1.005e+05 on 7 and 1103344 DF, p-value: < 2.2e-16 #with mat age category
m_BW1_matage_cat = update(m_BW1, . ~ . + mat_age_cat)
summary(m_BW1_matage_cat)##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex +
## mat_age_cat, data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2703.1 -264.2 -13.3 249.7 5611.5
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.714e+03 8.897e+00 -417.387 < 2e-16 ***
## GA_weeks 1.790e+02 2.241e-01 798.686 < 2e-16 ***
## birth_Y_M_num -1.055e-01 7.897e-03 -13.362 < 2e-16 ***
## parity_cat(1,2] 1.351e+02 8.333e-01 162.099 < 2e-16 ***
## parity_cat(2,3] 1.792e+02 1.275e+00 140.475 < 2e-16 ***
## parity_cat(3,20] 1.999e+02 2.217e+00 90.172 < 2e-16 ***
## sex2 -1.469e+02 7.515e-01 -195.549 < 2e-16 ***
## mat_age_cat(10,20] -2.391e+01 3.549e+00 -6.737 1.62e-11 ***
## mat_age_cat(20,25] 1.168e+00 1.401e+00 0.834 0.40440
## mat_age_cat(30,35] -3.325e+00 9.442e-01 -3.521 0.00043 ***
## mat_age_cat(35,40] -3.223e+00 1.118e+00 -2.883 0.00394 **
## mat_age_cat(40,70] -3.808e+00 2.071e+00 -1.839 0.06591 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 394.4 on 1103340 degrees of freedom
## Multiple R-squared: 0.3895, Adjusted R-squared: 0.3895
## F-statistic: 6.399e+04 on 11 and 1103340 DF, p-value: < 2.2e-16 #with mother nationality (categorical 1)
m_BW1_mat_nat_cat1 = update(m_BW1_matage_cat, . ~ . + mother_nationality_cat1)
summary(m_BW1_mat_nat_cat1)##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex +
## mat_age_cat + mother_nationality_cat1, data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2692.6 -263.9 -13.1 249.5 5599.7
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.726e+03 8.921e+00 -417.720 < 2e-16 ***
## GA_weeks 1.791e+02 2.246e-01 797.554 < 2e-16 ***
## birth_Y_M_num -1.158e-01 7.909e-03 -14.647 < 2e-16 ***
## parity_cat(1,2] 1.355e+02 8.342e-01 162.450 < 2e-16 ***
## parity_cat(2,3] 1.802e+02 1.277e+00 141.124 < 2e-16 ***
## parity_cat(3,20] 2.017e+02 2.221e+00 90.822 < 2e-16 ***
## sex2 -1.470e+02 7.524e-01 -195.377 < 2e-16 ***
## mat_age_cat(10,20] -2.610e+01 3.591e+00 -7.268 3.65e-13 ***
## mat_age_cat(20,25] -2.189e+00 1.411e+00 -1.551 0.1209
## mat_age_cat(30,35] -2.432e+00 9.457e-01 -2.572 0.0101 *
## mat_age_cat(35,40] -2.972e+00 1.119e+00 -2.656 0.0079 **
## mat_age_cat(40,70] -4.134e+00 2.072e+00 -1.995 0.0460 *
## mother_nationality_cat1Outside Switzerland 2.568e+01 7.791e-01 32.962 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 394.2 on 1099355 degrees of freedom
## (3984 observations effacées parce que manquantes)
## Multiple R-squared: 0.39, Adjusted R-squared: 0.39
## F-statistic: 5.857e+04 on 12 and 1099355 DF, p-value: < 2.2e-16 #with mother nationality (categorical 2)
m_BW1_mat_nat_cat2 = update(m_BW1_matage_cat, . ~ . + mother_nationality_cat2)
summary(m_BW1_mat_nat_cat2)##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex +
## mat_age_cat + mother_nationality_cat2, data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2692.9 -263.7 -13.2 249.4 5593.7
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.725e+03 8.921e+00 -417.531 < 2e-16 ***
## GA_weeks 1.790e+02 2.246e-01 797.302 < 2e-16 ***
## birth_Y_M_num -1.095e-01 7.912e-03 -13.843 < 2e-16 ***
## parity_cat(1,2] 1.359e+02 8.343e-01 162.851 < 2e-16 ***
## parity_cat(2,3] 1.813e+02 1.278e+00 141.884 < 2e-16 ***
## parity_cat(3,20] 2.040e+02 2.223e+00 91.768 < 2e-16 ***
## sex2 -1.470e+02 7.522e-01 -195.392 < 2e-16 ***
## mat_age_cat(10,20] -2.530e+01 3.591e+00 -7.045 1.86e-12 ***
## mat_age_cat(20,25] -1.942e+00 1.411e+00 -1.376 0.168890
## mat_age_cat(30,35] -2.970e+00 9.460e-01 -3.140 0.001692 **
## mat_age_cat(35,40] -3.931e+00 1.120e+00 -3.509 0.000449 ***
## mat_age_cat(40,70] -5.230e+00 2.072e+00 -2.524 0.011601 *
## mother_nationality_cat2Africa -4.329e+00 2.272e+00 -1.905 0.056735 .
## mother_nationality_cat2Asia -1.159e+01 2.067e+00 -5.608 2.05e-08 ***
## mother_nationality_cat2Europe 3.103e+01 8.486e-01 36.567 < 2e-16 ***
## mother_nationality_cat2Northern America 5.912e+01 4.885e+00 12.101 < 2e-16 ***
## mother_nationality_cat2Oceania 6.280e+01 1.122e+01 5.597 2.18e-08 ***
## mother_nationality_cat2Southern and Central America 4.453e+01 2.813e+00 15.832 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 394 on 1099310 degrees of freedom
## (4024 observations effacées parce que manquantes)
## Multiple R-squared: 0.3904, Adjusted R-squared: 0.3904
## F-statistic: 4.141e+04 on 17 and 1099310 DF, p-value: < 2.2e-16 #with urban category (categorical 3)
m_BW1_urb = update(m_BW1_mat_nat_cat2, . ~ . + Urban3)
summary(m_BW1_urb) ##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex +
## mat_age_cat + mother_nationality_cat2 + Urban3, data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2690.9 -263.7 -13.2 249.5 5596.4
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.726e+03 8.925e+00 -417.541 < 2e-16 ***
## GA_weeks 1.790e+02 2.246e-01 797.230 < 2e-16 ***
## birth_Y_M_num -1.098e-01 7.912e-03 -13.871 < 2e-16 ***
## parity_cat(1,2] 1.361e+02 8.354e-01 162.944 < 2e-16 ***
## parity_cat(2,3] 1.818e+02 1.280e+00 141.989 < 2e-16 ***
## parity_cat(3,20] 2.046e+02 2.225e+00 91.942 < 2e-16 ***
## sex2 -1.470e+02 7.522e-01 -195.388 < 2e-16 ***
## mat_age_cat(10,20] -2.537e+01 3.591e+00 -7.067 1.59e-12 ***
## mat_age_cat(20,25] -1.955e+00 1.411e+00 -1.385 0.165957
## mat_age_cat(30,35] -3.169e+00 9.466e-01 -3.348 0.000815 ***
## mat_age_cat(35,40] -4.288e+00 1.122e+00 -3.823 0.000132 ***
## mat_age_cat(40,70] -5.652e+00 2.073e+00 -2.726 0.006408 **
## mother_nationality_cat2Africa -5.606e+00 2.282e+00 -2.456 0.014036 *
## mother_nationality_cat2Asia -1.276e+01 2.076e+00 -6.145 8.00e-10 ***
## mother_nationality_cat2Europe 3.043e+01 8.545e-01 35.612 < 2e-16 ***
## mother_nationality_cat2Northern America 5.810e+01 4.888e+00 11.886 < 2e-16 ***
## mother_nationality_cat2Oceania 6.192e+01 1.122e+01 5.518 3.43e-08 ***
## mother_nationality_cat2Southern and Central America 4.354e+01 2.818e+00 15.452 < 2e-16 ***
## Urban3 4.557e+00 7.656e-01 5.952 2.64e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 394 on 1099309 degrees of freedom
## (4024 observations effacées parce que manquantes)
## Multiple R-squared: 0.3904, Adjusted R-squared: 0.3904
## F-statistic: 3.911e+04 on 18 and 1099309 DF, p-value: < 2.2e-16 #with mean ssep2
m_BW1_ssep = update(m_BW1_urb, . ~ . + mean_ssep2)
summary(m_BW1_ssep) ##
## Call:
## lm(formula = BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex +
## mat_age_cat + mother_nationality_cat2 + Urban3 + mean_ssep2,
## data = bevn_eco_in7)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2683.9 -263.6 -13.2 249.4 5637.2
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.812e+03 9.264e+00 -411.470 < 2e-16 ***
## GA_weeks 1.790e+02 2.245e-01 797.495 < 2e-16 ***
## birth_Y_M_num -1.056e-01 7.909e-03 -13.356 < 2e-16 ***
## parity_cat(1,2] 1.372e+02 8.355e-01 164.206 < 2e-16 ***
## parity_cat(2,3] 1.846e+02 1.282e+00 143.986 < 2e-16 ***
## parity_cat(3,20] 2.092e+02 2.228e+00 93.885 < 2e-16 ***
## sex2 -1.469e+02 7.518e-01 -195.401 < 2e-16 ***
## mat_age_cat(10,20] -2.287e+01 3.590e+00 -6.372 1.86e-10 ***
## mat_age_cat(20,25] 5.588e-02 1.412e+00 0.040 0.968
## mat_age_cat(30,35] -6.142e+00 9.500e-01 -6.465 1.01e-10 ***
## mat_age_cat(35,40] -9.222e+00 1.130e+00 -8.159 3.39e-16 ***
## mat_age_cat(40,70] -1.114e+01 2.078e+00 -5.358 8.41e-08 ***
## mother_nationality_cat2Africa -3.506e+00 2.282e+00 -1.536 0.124
## mother_nationality_cat2Asia -1.401e+01 2.075e+00 -6.751 1.47e-11 ***
## mother_nationality_cat2Europe 3.073e+01 8.541e-01 35.981 < 2e-16 ***
## mother_nationality_cat2Northern America 5.289e+01 4.888e+00 10.820 < 2e-16 ***
## mother_nationality_cat2Oceania 5.483e+01 1.122e+01 4.888 1.02e-06 ***
## mother_nationality_cat2Southern and Central America 4.300e+01 2.816e+00 15.269 < 2e-16 ***
## Urban3 -1.736e-01 7.776e-01 -0.223 0.823
## mean_ssep2 1.549e+00 4.528e-02 34.214 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 393.8 on 1099308 degrees of freedom
## (4024 observations effacées parce que manquantes)
## Multiple R-squared: 0.391, Adjusted R-squared: 0.391
## F-statistic: 3.715e+04 on 19 and 1099308 DF, p-value: < 2.2e-16models <- mtable("m1: GA + birthmonth + parity + sex"=m_BW1, "m2: m1 + mat.age cont."=m_BW1_matage_cont, "m3: m1 + mat.age.categ"=m_BW1_matage_cat, "m4: m3 + mother nationality_cat1"= m_BW1_mat_nat_cat1 , "m5: m3 + mother nationality_cat2"= m_BW1_mat_nat_cat2, "m6: m5+ Urban"=m_BW1_urb, "m7: m6+ ssep"=m_BW1_ssep, coef.style="ci.horizontal",
summary.stats=c("Estimate","Pr(>|t|)","AIC")
)
write_html(models, here("output", "summary_BW_lm_models.html"))
anova(m_BW1_mat_nat_cat2,m_BW1_urb)## Analysis of Variance Table
##
## Model 1: BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex + mat_age_cat +
## mother_nationality_cat2
## Model 2: BW ~ GA_weeks + birth_Y_M_num + parity_cat + sex + mat_age_cat +
## mother_nationality_cat2 + Urban3
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1099310 1.7069e+11
## 2 1099309 1.7068e+11 1 5501217 35.431 2.643e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Based on the last model (m_BW1_urb_cat_3): - GA: each gestational age increase of one week increase BW by 179g. - Birthdate (month + year) seems to have a significant but irrelevant effect (-0.11g each month) - Parity: with each category, the birthweight increases compared to primiparous women: +137g for parity=2, +185g for parity=3, +209g for parity >3. - sex: females are 147g lighter than males - maternal age (cat): compared to the reference category of 25-30yo, women of any other age category have lighter babies (all significant except for category 20-25), but the effect is worth noticing mainly for younger women (10-20), who have babies 23g lighter. - mother nationality: compared to babies born from mother who have the Swiss nationality, babies born from Asian mothers are a bit lighter (-14g), while babies born from Europe (excl Switz), Northern America, Oceania or Southern/Central America are bigger (+31g, +53, +55, +43g respectively). - Urbanicity: once SSEP is added to the model, being born in an urban area no longer increases by 4.5g birthweight.
Stillbirth
#baseline model (with maternal nationality categ1)
m_SB_multiv_1 <- glm(stillbirth~mat_age_cat + BW + GA_weeks + birth_Y_M_num + mother_nationality_cat1,
family= "binomial",data=bevn_eco_in6)
summary(m_SB_multiv_1) ##
## Call:
## glm(formula = stillbirth ~ mat_age_cat + BW + GA_weeks + birth_Y_M_num +
## mother_nationality_cat1, family = "binomial", data = bevn_eco_in6)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.2716 -0.0562 -0.0422 -0.0319 4.5593
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.056e+00 2.211e-01 18.346 < 2e-16 ***
## mat_age_cat(10,20] 2.232e-02 1.527e-01 0.146 0.88380
## mat_age_cat(20,25] 4.084e-02 7.015e-02 0.582 0.56042
## mat_age_cat(30,35] -6.876e-02 4.948e-02 -1.390 0.16460
## mat_age_cat(35,40] 1.429e-01 5.372e-02 2.660 0.00781 **
## mat_age_cat(40,70] 1.642e-01 8.746e-02 1.878 0.06037 .
## BW -1.633e-03 5.225e-05 -31.250 < 2e-16 ***
## GA_weeks -1.419e-01 9.291e-03 -15.273 < 2e-16 ***
## birth_Y_M_num -1.529e-04 4.001e-04 -0.382 0.70240
## mother_nationality_cat1Outside Switzerland 5.137e-02 3.873e-02 1.326 0.18473
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 44819 on 1102656 degrees of freedom
## Residual deviance: 30018 on 1102647 degrees of freedom
## (4022 observations effacées parce que manquantes)
## AIC: 30038
##
## Number of Fisher Scoring iterations: 9 exp(coef(summary(m_SB_multiv_1)))[,"Estimate"]## (Intercept) mat_age_cat(10,20]
## 57.7267534 1.0225664
## mat_age_cat(20,25] mat_age_cat(30,35]
## 1.0416895 0.9335470
## mat_age_cat(35,40] mat_age_cat(40,70]
## 1.1536083 1.1785077
## BW GA_weeks
## 0.9983684 0.8677084
## birth_Y_M_num mother_nationality_cat1Outside Switzerland
## 0.9998471 1.0527138#With mother nationality category 2
m_SB_multiv_2 = update(m_SB_multiv_1, . ~ . - mother_nationality_cat1 + mother_nationality_cat2)
summary(m_SB_multiv_2) ##
## Call:
## glm(formula = stillbirth ~ mat_age_cat + BW + GA_weeks + birth_Y_M_num +
## mother_nationality_cat2, family = "binomial", data = bevn_eco_in6)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3273 -0.0562 -0.0421 -0.0319 4.5580
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.058e+00 2.211e-01 18.353 < 2e-16 ***
## mat_age_cat(10,20] 1.988e-02 1.527e-01 0.130 0.89640
## mat_age_cat(20,25] 3.681e-02 7.018e-02 0.524 0.59998
## mat_age_cat(30,35] -6.613e-02 4.950e-02 -1.336 0.18154
## mat_age_cat(35,40] 1.481e-01 5.378e-02 2.754 0.00589 **
## mat_age_cat(40,70] 1.672e-01 8.753e-02 1.910 0.05613 .
## BW -1.631e-03 5.227e-05 -31.197 < 2e-16 ***
## GA_weeks -1.421e-01 9.292e-03 -15.294 < 2e-16 ***
## birth_Y_M_num -1.689e-04 4.003e-04 -0.422 0.67313
## mother_nationality_cat2Africa 1.968e-01 9.332e-02 2.109 0.03498 *
## mother_nationality_cat2Asia -5.127e-02 1.015e-01 -0.505 0.61339
## mother_nationality_cat2Europe 6.090e-02 4.254e-02 1.432 0.15225
## mother_nationality_cat2Northern America -3.749e-01 2.891e-01 -1.297 0.19465
## mother_nationality_cat2Oceania -1.405e-01 7.324e-01 -0.192 0.84790
## mother_nationality_cat2Southern and Central America -7.066e-02 1.361e-01 -0.519 0.60375
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 44818 on 1102616 degrees of freedom
## Residual deviance: 30010 on 1102602 degrees of freedom
## (4062 observations effacées parce que manquantes)
## AIC: 30040
##
## Number of Fisher Scoring iterations: 9 exp(coef(summary(m_SB_multiv_2)))[,"Estimate"]## (Intercept) mat_age_cat(10,20]
## 57.8337897 1.0200823
## mat_age_cat(20,25] mat_age_cat(30,35]
## 1.0374924 0.9360100
## mat_age_cat(35,40] mat_age_cat(40,70]
## 1.1596379 1.1819638
## BW GA_weeks
## 0.9983707 0.8675250
## birth_Y_M_num mother_nationality_cat2Africa
## 0.9998311 1.2174615
## mother_nationality_cat2Asia mother_nationality_cat2Europe
## 0.9500227 1.0627909
## mother_nationality_cat2Northern America mother_nationality_cat2Oceania
## 0.6873574 0.8689497
## mother_nationality_cat2Southern and Central America
## 0.9317790#adding Urbanicity, language and ssep
m_SB_multiv_3 = update(m_SB_multiv_2, . ~ . + Language + mean_ssep2)
summary(m_SB_multiv_3) ##
## Call:
## glm(formula = stillbirth ~ mat_age_cat + BW + GA_weeks + birth_Y_M_num +
## mother_nationality_cat2 + Language + mean_ssep2, family = "binomial",
## data = bevn_eco_in6)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3414 -0.0561 -0.0421 -0.0319 4.5727
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.356e+00 2.605e-01 16.721 < 2e-16 ***
## mat_age_cat(10,20] 1.722e-02 1.527e-01 0.113 0.91021
## mat_age_cat(20,25] 3.171e-02 7.021e-02 0.452 0.65152
## mat_age_cat(30,35] -5.716e-02 4.964e-02 -1.152 0.24951
## mat_age_cat(35,40] 1.622e-01 5.413e-02 2.996 0.00273 **
## mat_age_cat(40,70] 1.820e-01 8.786e-02 2.071 0.03837 *
## BW -1.625e-03 5.233e-05 -31.046 < 2e-16 ***
## GA_weeks -1.430e-01 9.303e-03 -15.370 < 2e-16 ***
## birth_Y_M_num -1.975e-04 4.007e-04 -0.493 0.62206
## mother_nationality_cat2Africa 1.817e-01 9.408e-02 1.931 0.05345 .
## mother_nationality_cat2Asia -4.477e-02 1.016e-01 -0.441 0.65938
## mother_nationality_cat2Europe 6.008e-02 4.261e-02 1.410 0.15849
## mother_nationality_cat2Northern America -3.626e-01 2.895e-01 -1.252 0.21041
## mother_nationality_cat2Oceania -1.144e-01 7.317e-01 -0.156 0.87575
## mother_nationality_cat2Southern and Central America -6.794e-02 1.364e-01 -0.498 0.61840
## LanguageFrench 2.668e-02 4.435e-02 0.602 0.54739
## LanguageItalian -8.660e-02 1.027e-01 -0.843 0.39909
## mean_ssep2 -4.997e-03 2.308e-03 -2.166 0.03035 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 44818 on 1102616 degrees of freedom
## Residual deviance: 30004 on 1102599 degrees of freedom
## (4062 observations effacées parce que manquantes)
## AIC: 30040
##
## Number of Fisher Scoring iterations: 9 exp(coef(summary(m_SB_multiv_3)))[,"Estimate"]## (Intercept) mat_age_cat(10,20]
## 77.9099925 1.0173651
## mat_age_cat(20,25] mat_age_cat(30,35]
## 1.0322182 0.9444393
## mat_age_cat(35,40] mat_age_cat(40,70]
## 1.1760913 1.1995569
## BW GA_weeks
## 0.9983766 0.8667630
## birth_Y_M_num mother_nationality_cat2Africa
## 0.9998025 1.1992535
## mother_nationality_cat2Asia mother_nationality_cat2Europe
## 0.9562133 1.0619240
## mother_nationality_cat2Northern America mother_nationality_cat2Oceania
## 0.6958439 0.8918948
## mother_nationality_cat2Southern and Central America LanguageFrench
## 0.9343119 1.0270414
## LanguageItalian mean_ssep2
## 0.9170481 0.9950154 exp(confint(m_SB_multiv_3))## Attente de la réalisation du profilage...## 2.5 % 97.5 %
## (Intercept) 46.7927092 129.9148522
## mat_age_cat(10,20] 0.7469331 1.3600562
## mat_age_cat(20,25] 0.8985277 1.1833121
## mat_age_cat(30,35] 0.8569724 1.0411078
## mat_age_cat(35,40] 1.0575927 1.3076687
## mat_age_cat(40,70] 1.0073739 1.4218207
## BW 0.9982742 0.9984790
## GA_weeks 0.8510620 0.8826742
## birth_Y_M_num 0.9990177 1.0005883
## mother_nationality_cat2Africa 0.9940257 1.4376582
## mother_nationality_cat2Asia 0.7799863 1.1618574
## mother_nationality_cat2Europe 0.9765468 1.1540946
## mother_nationality_cat2Northern America 0.3772281 1.1823654
## mother_nationality_cat2Oceania 0.1434846 2.9190497
## mother_nationality_cat2Southern and Central America 0.7090464 1.2111531
## LanguageFrench 0.9412045 1.1199380
## LanguageItalian 0.7462825 1.1165289
## mean_ssep2 0.9905240 0.9995256Based on the last model (m_SB_multiv_3): - maternal age: compared to age 25-30, being over 35 increases the risk of having a stillborn baby: CI95% 1.06-1.31 for mothers 35-40yo, and CI95% 1.01-1.42 for mothers 40-70yo. - increasing BW slightly decreases stillborn risk: CI95% 0.9983-0.9985 - Gestational age: increasing GA decreases risk of stillborn: CI95% 0.85-0.88 - Mother nationality: compared to babies born from Swiss-mothers, babies born from mothers with African nationality have an increased risk of being stillborn (pvalue borderline = 0.053), CI95% 0.99-1.44 - SSEP of the place of residence of the mother plays a minor role: the higher the ssep, the lower the risk of having a stillborn baby: CI95% 0.990-0.999 per ssep point increase
GAM Models
Birthweight
2007-2020
m_BW_bam4 <- bam(BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat1 + Language, data=bevn_eco_in7)
summary(m_BW_bam4)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat1 +
## Language
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3338.3982 0.8351 3997.771 <2e-16 ***
## parity_cat(1,2] 127.1812 0.8344 152.424 <2e-16 ***
## parity_cat(2,3] 172.7037 1.2761 135.337 <2e-16 ***
## parity_cat(3,20] 197.8559 2.2120 89.445 <2e-16 ***
## sex2 -147.0453 0.7460 -197.117 <2e-16 ***
## Urban3 -0.1677 0.7896 -0.212 0.832
## mother_nationality_cat1Outside Switzerland 27.7466 0.7860 35.301 <2e-16 ***
## LanguageFrench -47.8511 0.9157 -52.259 <2e-16 ***
## LanguageItalian -62.3288 2.0496 -30.410 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.865 8.993 73229.115 < 2e-16 ***
## s(mat_age) 5.982 6.893 14.740 < 2e-16 ***
## s(birth_Y_M_num) 8.276 8.856 40.858 < 2e-16 ***
## s(month) 5.905 7.065 3.169 0.00228 **
## s(mean_ssep2) 7.730 8.592 44.345 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.4 Deviance explained = 40%
## fREML = 8.1211e+06 Scale est. = 1.5269e+05 n = 1099368# with mat. nationality category 2
m_BW_bam5 <- bam(BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7)
summary(m_BW_bam5)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3335.2692 0.8361 3988.893 < 2e-16 ***
## parity_cat(1,2] 127.4690 0.8345 152.757 < 2e-16 ***
## parity_cat(2,3] 173.6228 1.2772 135.937 < 2e-16 ***
## parity_cat(3,20] 199.6926 2.2149 90.158 < 2e-16 ***
## sex2 -147.0105 0.7458 -197.129 < 2e-16 ***
## Urban3 0.6597 0.7906 0.834 0.404
## LanguageFrench -48.0928 0.9182 -52.375 < 2e-16 ***
## LanguageItalian -63.1908 2.0493 -30.835 < 2e-16 ***
## mother_nationality_cat2Africa 11.0218 2.2735 4.848 1.25e-06 ***
## mother_nationality_cat2Asia -17.9221 2.0596 -8.702 < 2e-16 ***
## mother_nationality_cat2Europe 32.5723 0.8501 38.315 < 2e-16 ***
## mother_nationality_cat2Northern America 65.0137 4.8546 13.392 < 2e-16 ***
## mother_nationality_cat2Oceania 63.4769 11.1288 5.704 1.17e-08 ***
## mother_nationality_cat2Southern and Central America 47.0430 2.7991 16.806 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.825 8.987 73218.407 <2e-16 ***
## s(mat_age) 5.994 6.904 16.693 <2e-16 ***
## s(birth_Y_M_num) 8.285 8.859 39.007 <2e-16 ***
## s(month) 5.895 7.054 3.183 0.0022 **
## s(mean_ssep2) 7.634 8.534 43.681 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.401 Deviance explained = 40.1%
## fREML = 8.1205e+06 Scale est. = 1.5259e+05 n = 1099328 plot(m_BW_bam5, select=1, ylim=c(0, 4000),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam5)[1])
plot(m_BW_bam5, select=2, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam5)[1])
plot(m_BW_bam5, select=3, ylim=c(3300, 3380),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam5)[1])
plot(m_BW_bam5, select=4, ylim=c(3325, 3345),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam5)[1])
plot(m_BW_bam5, select=5, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam5)[1])
plot(m_BW_bam5, all.terms = TRUE,pages=1, shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam5)[1])
#checking that residuals are not determined
gam.check(m_BW_bam5)
##
## Method: fREML Optimizer: perf newton
## full convergence after 15 iterations.
## Gradient range [-0.09472596,0.06729881]
## (score 8120451 & scale 152594.9).
## Hessian positive definite, eigenvalue range [1.614702,549654.6].
## Model rank = 59 / 59
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.83 1.00 0.45
## s(mat_age) 9.00 5.99 0.99 0.15
## s(birth_Y_M_num) 9.00 8.28 0.99 0.21
## s(month) 9.00 5.89 0.99 0.14
## s(mean_ssep2) 9.00 7.63 1.00 0.58 concurvity(m_BW_bam5, full=TRUE) ## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.8015821 0.017148348 0.12307078 0.030160407 0.028069842 0.17908052
## observed 0.8015821 0.006380490 0.04229698 0.007225289 0.002128965 0.10239571
## estimate 0.8015821 0.009738434 0.08228622 0.011467296 0.021047680 0.08954061 #Interpretation: p-values of the GAM check are not significant: residuals are randomly distributed. But k-index is <1 for mat age and month
#try with maternal age with k=12
m_BW_bam5_1 <- bam(BW ~ s(GA_weeks) + s(mat_age, k=13) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7)
summary(m_BW_bam5_1)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age, k = 13) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3335.2688 0.8361 3988.887 < 2e-16 ***
## parity_cat(1,2] 127.4691 0.8345 152.758 < 2e-16 ***
## parity_cat(2,3] 173.6230 1.2772 135.937 < 2e-16 ***
## parity_cat(3,20] 199.6909 2.2149 90.157 < 2e-16 ***
## sex2 -147.0105 0.7458 -197.129 < 2e-16 ***
## Urban3 0.6594 0.7906 0.834 0.404
## LanguageFrench -48.0927 0.9182 -52.374 < 2e-16 ***
## LanguageItalian -63.1902 2.0493 -30.835 < 2e-16 ***
## mother_nationality_cat2Africa 11.0235 2.2735 4.849 1.24e-06 ***
## mother_nationality_cat2Asia -17.9212 2.0596 -8.701 < 2e-16 ***
## mother_nationality_cat2Europe 32.5733 0.8501 38.316 < 2e-16 ***
## mother_nationality_cat2Northern America 65.0142 4.8546 13.392 < 2e-16 ***
## mother_nationality_cat2Oceania 63.4731 11.1288 5.704 1.17e-08 ***
## mother_nationality_cat2Southern and Central America 47.0407 2.7991 16.806 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.828 8.988 73215.221 <2e-16 ***
## s(mat_age) 6.568 7.836 14.641 <2e-16 ***
## s(birth_Y_M_num) 8.285 8.859 39.006 <2e-16 ***
## s(month) 5.893 7.053 3.183 0.0022 **
## s(mean_ssep2) 7.634 8.534 43.676 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.401 Deviance explained = 40.1%
## fREML = 8.1205e+06 Scale est. = 1.5259e+05 n = 1099328 gam.check(m_BW_bam5_1) #pvalue not significant, k-index for month and birthdate is < 1.
##
## Method: fREML Optimizer: perf newton
## full convergence after 12 iterations.
## Gradient range [-6.281771e-07,8.370701e-08]
## (score 8120452 & scale 152594.9).
## Hessian positive definite, eigenvalue range [1.614311,549654.5].
## Model rank = 62 / 62
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.83 1.01 0.66
## s(mat_age) 12.00 6.57 1.01 0.71
## s(birth_Y_M_num) 9.00 8.28 0.99 0.23
## s(month) 9.00 5.89 1.00 0.54
## s(mean_ssep2) 9.00 7.63 1.01 0.68#try with mat age k=12 and month k=10
m_BW_bam5_2 <- bam(BW ~ s(GA_weeks) + s(mat_age, k=13) + s(birth_Y_M_num) + s(month, k=11) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7)
summary(m_BW_bam5_2)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age, k = 13) + s(birth_Y_M_num) + s(month,
## k = 11) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language +
## mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3335.2684 0.8361 3988.889 < 2e-16 ***
## parity_cat(1,2] 127.4692 0.8345 152.758 < 2e-16 ***
## parity_cat(2,3] 173.6232 1.2772 135.937 < 2e-16 ***
## parity_cat(3,20] 199.6911 2.2149 90.157 < 2e-16 ***
## sex2 -147.0104 0.7458 -197.129 < 2e-16 ***
## Urban3 0.6595 0.7906 0.834 0.404
## LanguageFrench -48.0926 0.9182 -52.374 < 2e-16 ***
## LanguageItalian -63.1901 2.0493 -30.835 < 2e-16 ***
## mother_nationality_cat2Africa 11.0235 2.2735 4.849 1.24e-06 ***
## mother_nationality_cat2Asia -17.9212 2.0596 -8.702 < 2e-16 ***
## mother_nationality_cat2Europe 32.5734 0.8501 38.316 < 2e-16 ***
## mother_nationality_cat2Northern America 65.0139 4.8546 13.392 < 2e-16 ***
## mother_nationality_cat2Oceania 63.4734 11.1288 5.704 1.17e-08 ***
## mother_nationality_cat2Southern and Central America 47.0408 2.7991 16.806 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.828 8.988 73215.198 < 2e-16 ***
## s(mat_age) 6.568 7.836 14.641 < 2e-16 ***
## s(birth_Y_M_num) 8.285 8.859 39.006 < 2e-16 ***
## s(month) 5.978 7.248 3.072 0.00277 **
## s(mean_ssep2) 7.634 8.534 43.676 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.401 Deviance explained = 40.1%
## fREML = 8.1205e+06 Scale est. = 1.5259e+05 n = 1099328 gam.check(m_BW_bam5_2) 
##
## Method: fREML Optimizer: perf newton
## full convergence after 11 iterations.
## Gradient range [-0.0009357942,0.0007208003]
## (score 8120452 & scale 152595).
## Hessian positive definite, eigenvalue range [1.648423,549654.5].
## Model rank = 63 / 63
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.83 0.99 0.17
## s(mat_age) 12.00 6.57 0.98 0.11
## s(birth_Y_M_num) 9.00 8.28 1.00 0.50
## s(month) 10.00 5.98 0.99 0.20
## s(mean_ssep2) 9.00 7.63 1.00 0.52concurvity(m_BW_bam5_2, full=TRUE) ## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.8015963 0.017211066 0.12307439 0.030165524 0.028070793 0.17908212
## observed 0.8015963 0.006386136 0.04255397 0.007228201 0.002163185 0.10239354
## estimate 0.8015963 0.009760255 0.07758667 0.011469079 0.019775195 0.08954325 #from the "worst" case, the values are not high (high would be >0.8, usually)
plot(m_BW_bam5_2, select=3, ylim=c(-15, 20), xlim=c(144, 168),shade = TRUE, shade.col = "lightblue")
#comparing models with different k parameters
AIC(m_BW_bam5, m_BW_bam5_1)## df AIC
## m_BW_bam5 53.45272 16240889
## m_BW_bam5_1 54.04528 16240890 AIC(m_BW_bam5_2, m_BW_bam5_1)## df AIC
## m_BW_bam5_2 54.13589 16240891
## m_BW_bam5_1 54.04528 16240890 #with mat nationality category 3
m_BW_bam6 <- bam(BW ~ s(GA_weeks) + s(mat_age, k=13) + s(birth_Y_M_num) + s(month, k=11) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat3, data=bevn_eco_in7)
summary(m_BW_bam6)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age, k = 13) + s(birth_Y_M_num) + s(month,
## k = 11) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language +
## mother_nationality_cat3
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3336.7327 0.8368 3987.583 <2e-16 ***
## parity_cat(1,2] 125.8337 0.8356 150.600 <2e-16 ***
## parity_cat(2,3] 169.0143 1.2812 131.915 <2e-16 ***
## parity_cat(3,20] 193.4530 2.2145 87.357 <2e-16 ***
## sex2 -146.9714 0.7446 -197.379 <2e-16 ***
## Urban3 0.1904 0.7896 0.241 0.809
## LanguageFrench -38.6042 0.9461 -40.805 <2e-16 ***
## LanguageItalian -56.7421 2.0605 -27.538 <2e-16 ***
## mother_nationality_cat3Central Europe 44.5600 1.5104 29.503 <2e-16 ***
## mother_nationality_cat3Eastern Europe 102.0017 3.7786 26.995 <2e-16 ***
## mother_nationality_cat3Northern Europe 89.5490 5.9617 15.021 <2e-16 ***
## mother_nationality_cat3other 12.4683 1.3318 9.362 <2e-16 ***
## mother_nationality_cat3SE Europe 73.8194 1.4196 51.999 <2e-16 ***
## mother_nationality_cat3Southern Europe 2.4753 1.8727 1.322 0.186
## mother_nationality_cat3SW Europe -46.9685 1.8041 -26.034 <2e-16 ***
## mother_nationality_cat3Western Europe 24.1616 2.1979 10.993 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.866 8.993 73364.658 <2e-16 ***
## s(mat_age) 5.994 7.251 11.542 <2e-16 ***
## s(birth_Y_M_num) 8.311 8.868 40.084 <2e-16 ***
## s(month) 5.937 7.203 3.075 0.0028 **
## s(mean_ssep2) 7.410 8.390 36.484 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.403 Deviance explained = 40.3%
## fREML = 8.1191e+06 Scale est. = 1.5213e+05 n = 10993682008-2010
m_BW_bam_2009 <- bam(BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7_2008_10)
summary(m_BW_bam_2009)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3329.5320 1.8664 1783.903 < 2e-16 ***
## parity_cat(1,2] 122.2021 1.8849 64.834 < 2e-16 ***
## parity_cat(2,3] 172.3539 2.8839 59.763 < 2e-16 ***
## parity_cat(3,20] 196.5353 4.9231 39.921 < 2e-16 ***
## sex2 -146.2370 1.6800 -87.046 < 2e-16 ***
## Urban3 0.8105 1.7448 0.465 0.64227
## LanguageFrench -46.0005 2.0609 -22.320 < 2e-16 ***
## LanguageItalian -65.1682 4.3466 -14.993 < 2e-16 ***
## mother_nationality_cat2Africa 41.2553 6.0376 6.833 8.33e-12 ***
## mother_nationality_cat2Asia -15.2957 4.9412 -3.096 0.00196 **
## mother_nationality_cat2Europe 36.4595 1.9397 18.797 < 2e-16 ***
## mother_nationality_cat2Northern America 82.8381 10.6237 7.797 6.34e-15 ***
## mother_nationality_cat2Oceania 54.2201 24.9694 2.171 0.02990 *
## mother_nationality_cat2Southern and Central America 60.6176 6.0057 10.093 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.245 8.807 14715.322 < 2e-16 ***
## s(mat_age) 3.355 4.179 6.894 9.88e-06 ***
## s(birth_Y_M_num) 2.494 3.063 17.617 < 2e-16 ***
## s(month) 2.167 2.695 1.748 0.142
## s(mean_ssep2) 1.005 1.010 70.227 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.397 Deviance explained = 39.7%
## fREML = 1.6251e+06 Scale est. = 1.548e+05 n = 219790 plot(m_BW_bam_2009, select=1, ylim=c(0, 4000),shade = TRUE, shade.col = "lightblue",shift = coef(m_BW_bam_2009)[1])
plot(m_BW_bam_2009, select=2, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue",shift = coef(m_BW_bam_2009)[1])
plot(m_BW_bam_2009, select=3, ylim=c(3310, 3360),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_2009)[1])
plot(m_BW_bam_2009, select=4, ylim=c(3310, 3340),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_2009)[1])
plot(m_BW_bam_2009, select=5, ylim=c(3300, 3400),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_2009)[1])
gam.check(m_BW_bam_2009)
##
## Method: fREML Optimizer: perf newton
## full convergence after 9 iterations.
## Gradient range [-0.005345534,0.002813906]
## (score 1625078 & scale 154799.6).
## Hessian positive definite, eigenvalue range [0.001879477,109885.5].
## Model rank = 59 / 59
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.24 1.01 0.73
## s(mat_age) 9.00 3.35 0.99 0.17
## s(birth_Y_M_num) 9.00 2.49 1.00 0.68
## s(month) 9.00 2.17 0.99 0.12
## s(mean_ssep2) 9.00 1.01 1.00 0.49 concurvity(m_BW_bam_2009, full=TRUE) #problem with birth_Y_M_num and month variables: values >>0.8, they correlate too much## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.7986046 0.016120856 0.13818504 0.9986011 0.9986011 0.18028011
## observed 0.7986046 0.006590663 0.02794266 0.1004240 0.8761449 0.13962224
## estimate 0.7986046 0.009381966 0.08881877 0.0985319 0.7612626 0.08259965 #removing month variable
m_BW_bam_2009_2 <- bam(BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7_2008_10)
gam.check(m_BW_bam_2009_2)
##
## Method: fREML Optimizer: perf newton
## full convergence after 9 iterations.
## Gradient range [-0.006332637,0.003326314]
## (score 1625081 & scale 154802.2).
## Hessian positive definite, eigenvalue range [0.00189899,109886].
## Model rank = 50 / 50
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.25 0.98 0.15
## s(mat_age) 9.00 3.35 0.99 0.15
## s(birth_Y_M_num) 9.00 2.84 1.01 0.70
## s(mean_ssep2) 9.00 1.01 1.00 0.57 concurvity(m_BW_bam_2009_2, full=TRUE) #values are now fine.## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(mean_ssep2)
## worst 0.798497 0.016067253 0.13813952 0.005952865 0.18024267
## observed 0.798497 0.006496112 0.02809684 0.005531003 0.13959235
## estimate 0.798497 0.009299257 0.08876812 0.003408024 0.08256354 summary(m_BW_bam_2009_2)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(mean_ssep2) +
## parity_cat + sex + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3329.564 1.866 1784.047 < 2e-16 ***
## parity_cat(1,2] 122.153 1.885 64.812 < 2e-16 ***
## parity_cat(2,3] 172.293 2.884 59.745 < 2e-16 ***
## parity_cat(3,20] 196.559 4.923 39.926 < 2e-16 ***
## sex2 -146.221 1.680 -87.037 < 2e-16 ***
## Urban3 0.818 1.745 0.469 0.63921
## LanguageFrench -45.975 2.061 -22.309 < 2e-16 ***
## LanguageItalian -65.095 4.346 -14.976 < 2e-16 ***
## mother_nationality_cat2Africa 41.267 6.038 6.835 8.23e-12 ***
## mother_nationality_cat2Asia -15.251 4.941 -3.087 0.00203 **
## mother_nationality_cat2Europe 36.459 1.940 18.797 < 2e-16 ***
## mother_nationality_cat2Northern America 82.856 10.624 7.799 6.26e-15 ***
## mother_nationality_cat2Oceania 54.309 24.970 2.175 0.02963 *
## mother_nationality_cat2Southern and Central America 60.686 6.006 10.105 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.246 8.808 14714.10 < 2e-16 ***
## s(mat_age) 3.354 4.178 6.89 9.95e-06 ***
## s(birth_Y_M_num) 2.838 3.533 14.97 < 2e-16 ***
## s(mean_ssep2) 1.005 1.010 70.18 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.397 Deviance explained = 39.7%
## fREML = 1.6251e+06 Scale est. = 1.548e+05 n = 219790 plot(m_BW_bam_2009_2, select=3, ylim=c(3310, 3360),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_2009_2)[1])
2017-2020
m_BW_bam_covid <- bam(BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7_2017_20)
summary(m_BW_bam_covid)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3341.837 1.513 2209.308 < 2e-16 ***
## parity_cat(1,2] 132.567 1.504 88.170 < 2e-16 ***
## parity_cat(2,3] 180.801 2.293 78.844 < 2e-16 ***
## parity_cat(3,20] 202.249 4.002 50.534 < 2e-16 ***
## sex2 -144.735 1.347 -107.482 < 2e-16 ***
## Urban3 -0.109 1.423 -0.077 0.9389
## LanguageFrench -51.202 1.633 -31.357 < 2e-16 ***
## LanguageItalian -62.481 3.949 -15.824 < 2e-16 ***
## mother_nationality_cat2Africa -6.570 3.768 -1.744 0.0812 .
## mother_nationality_cat2Asia -16.432 3.546 -4.634 3.59e-06 ***
## mother_nationality_cat2Europe 27.175 1.522 17.852 < 2e-16 ***
## mother_nationality_cat2Northern America 59.906 9.028 6.636 3.24e-11 ***
## mother_nationality_cat2Oceania 60.999 22.286 2.737 0.0062 **
## mother_nationality_cat2Southern and Central America 38.194 5.367 7.117 1.11e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.501 8.905 22220.130 <2e-16 ***
## s(mat_age) 5.972 6.889 19.738 <2e-16 ***
## s(birth_Y_M_num) 5.407 6.274 9.367 <2e-16 ***
## s(month) 4.886 5.963 2.472 0.0228 *
## s(mean_ssep2) 5.065 6.233 14.994 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.399 Deviance explained = 39.9%
## fREML = 2.4557e+06 Scale est. = 1.5056e+05 n = 332752 plot(m_BW_bam_covid, select=1, ylim=c(0, 4000),shade = TRUE, shade.col = "lightblue",shift = coef(m_BW_bam_covid)[1])
plot(m_BW_bam_covid, select=2, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue",shift = coef(m_BW_bam_covid)[1])
plot(m_BW_bam_covid, select=3, ylim=c(3310, 3360),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_covid)[1])
plot(m_BW_bam_covid, select=4, ylim=c(3330, 3360),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_covid)[1])
plot(m_BW_bam_covid, select=5, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_covid)[1])
gam.check(m_BW_bam_covid)
##
## Method: fREML Optimizer: perf newton
## full convergence after 9 iterations.
## Gradient range [-0.02973805,0.01546078]
## (score 2455711 & scale 150560).
## Hessian positive definite, eigenvalue range [0.8451653,166366.5].
## Model rank = 59 / 59
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.50 1.03 0.97
## s(mat_age) 9.00 5.97 1.00 0.51
## s(birth_Y_M_num) 9.00 5.41 1.02 0.88
## s(month) 9.00 4.89 1.02 0.85
## s(mean_ssep2) 9.00 5.06 1.01 0.67 concurvity(m_BW_bam_covid, full=TRUE)## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.8025327 0.018835869 0.10877183 0.98667852 0.98667868 0.17782830
## observed 0.8025327 0.007301885 0.06036586 0.04440909 0.08610047 0.08806201
## estimate 0.8025327 0.010995269 0.06973109 0.04756649 0.71521941 0.08027400 #removing month variable
m_BW_bam_covid_2 <- bam(BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7_2017_20)
summary(m_BW_bam_covid)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3341.837 1.513 2209.308 < 2e-16 ***
## parity_cat(1,2] 132.567 1.504 88.170 < 2e-16 ***
## parity_cat(2,3] 180.801 2.293 78.844 < 2e-16 ***
## parity_cat(3,20] 202.249 4.002 50.534 < 2e-16 ***
## sex2 -144.735 1.347 -107.482 < 2e-16 ***
## Urban3 -0.109 1.423 -0.077 0.9389
## LanguageFrench -51.202 1.633 -31.357 < 2e-16 ***
## LanguageItalian -62.481 3.949 -15.824 < 2e-16 ***
## mother_nationality_cat2Africa -6.570 3.768 -1.744 0.0812 .
## mother_nationality_cat2Asia -16.432 3.546 -4.634 3.59e-06 ***
## mother_nationality_cat2Europe 27.175 1.522 17.852 < 2e-16 ***
## mother_nationality_cat2Northern America 59.906 9.028 6.636 3.24e-11 ***
## mother_nationality_cat2Oceania 60.999 22.286 2.737 0.0062 **
## mother_nationality_cat2Southern and Central America 38.194 5.367 7.117 1.11e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.501 8.905 22220.130 <2e-16 ***
## s(mat_age) 5.972 6.889 19.738 <2e-16 ***
## s(birth_Y_M_num) 5.407 6.274 9.367 <2e-16 ***
## s(month) 4.886 5.963 2.472 0.0228 *
## s(mean_ssep2) 5.065 6.233 14.994 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.399 Deviance explained = 39.9%
## fREML = 2.4557e+06 Scale est. = 1.5056e+05 n = 332752 plot(m_BW_bam_covid_2, select=1, ylim=c(0, 4000),shade = TRUE, shade.col = "lightblue",shift = coef(m_BW_bam_covid_2)[1])
plot(m_BW_bam_covid_2, select=2, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue",shift = coef(m_BW_bam_covid_2)[1])
plot(m_BW_bam_covid_2, select=3, ylim=c(3310, 3360),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_covid_2)[1])
plot(m_BW_bam_covid_2, select=4, ylim=c(3200, 3400),shade = TRUE, shade.col = "lightblue", shift = coef(m_BW_bam_covid_2)[1])
gam.check(m_BW_bam_covid_2)
##
## Method: fREML Optimizer: perf newton
## full convergence after 12 iterations.
## Gradient range [-0.02525633,0.01270513]
## (score 2455714 & scale 150566.6).
## Hessian positive definite, eigenvalue range [0.9639783,166367].
## Model rank = 50 / 50
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.50 0.97 0.04 *
## s(mat_age) 9.00 5.97 1.01 0.83
## s(birth_Y_M_num) 9.00 5.64 0.99 0.15
## s(mean_ssep2) 9.00 5.07 1.00 0.46
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 concurvity(m_BW_bam_covid_2, full=TRUE) #values are now fine.## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(mean_ssep2)
## worst 0.8025004 0.018718378 0.10868673 0.0012305866 0.17778115
## observed 0.8025004 0.007202263 0.06041256 0.0006853439 0.08805890
## estimate 0.8025004 0.010886632 0.06965186 0.0009594351 0.08022341Stillbirth as outcome variable
2007-2020
m_SB_bam<-bam(stillbirth~s(GA_weeks)+ s(mat_age) + s(birth_Y_M_num) + s(month) +s(mean_ssep2) + Urban3 + Language + mother_nationality_cat2,family="binomial", data=bevn_eco_in6)
summary(m_SB_bam) #birth_Y_M_num has no effect##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.96844 0.04200 -165.927 <2e-16 ***
## Urban3 0.03816 0.03892 0.980 0.3269
## LanguageFrench 0.08108 0.04404 1.841 0.0656 .
## LanguageItalian -0.02042 0.10239 -0.199 0.8420
## mother_nationality_cat2Africa 0.22796 0.09350 2.438 0.0148 *
## mother_nationality_cat2Asia -0.03662 0.10164 -0.360 0.7187
## mother_nationality_cat2Europe 0.03212 0.04252 0.755 0.4500
## mother_nationality_cat2Northern America -0.34845 0.28509 -1.222 0.2216
## mother_nationality_cat2Oceania -0.17323 0.72070 -0.240 0.8100
## mother_nationality_cat2Southern and Central America -0.09875 0.13512 -0.731 0.4649
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.589 8.143 16734.682 < 2e-16 ***
## s(mat_age) 3.227 4.052 13.034 0.01214 *
## s(birth_Y_M_num) 1.013 1.026 0.010 0.95230
## s(month) 2.627 3.267 7.214 0.07017 .
## s(mean_ssep2) 1.008 1.016 6.865 0.00925 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.133 Deviance explained = 32.1%
## fREML = 1.5655e+06 Scale est. = 1 n = 1102617 plot(m_SB_bam, pages = 1, trans = plogis)
#adding intercept
plot(m_SB_bam, pages = 1, trans = plogis, shift = coef(m_SB_bam)[1])
# Second model_ without date of birth (Only month), and putting mean_ssep2 as linear term
m_SB_bam2 = update(m_SB_bam, . ~ . - s(mean_ssep2) + mean_ssep2 - s(birth_Y_M_num))
summary(m_SB_bam2)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(month) + Urban3 + Language +
## mother_nationality_cat2 + mean_ssep2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.614331 0.141044 -46.896 < 2e-16 ***
## Urban3 0.038204 0.038917 0.982 0.32626
## LanguageFrench 0.081151 0.044036 1.843 0.06535 .
## LanguageItalian -0.020692 0.102355 -0.202 0.83979
## mother_nationality_cat2Africa 0.228259 0.093470 2.442 0.01460 *
## mother_nationality_cat2Asia -0.036356 0.101616 -0.358 0.72051
## mother_nationality_cat2Europe 0.032281 0.042502 0.760 0.44754
## mother_nationality_cat2Northern America -0.348623 0.285086 -1.223 0.22138
## mother_nationality_cat2Oceania -0.173444 0.720702 -0.241 0.80982
## mother_nationality_cat2Southern and Central America -0.098668 0.135115 -0.730 0.46523
## mean_ssep2 -0.006102 0.002331 -2.617 0.00886 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.589 8.143 16736.825 <2e-16 ***
## s(mat_age) 3.224 4.049 13.026 0.0121 *
## s(month) 2.627 3.268 7.209 0.0703 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.133 Deviance explained = 32.1%
## fREML = 1.5655e+06 Scale est. = 1 n = 1102617 plot(m_SB_bam2, pages = 1, trans = plogis)
#adding intercept
plot(m_SB_bam2, pages = 1, trans = plogis, shift = coef(m_SB_bam)[1])
plot(m_SB_bam2, trans = plogis, select=1, shift = coef(m_SB_bam2)[1], shade = TRUE, shade.col = "lightblue") #GA
plot(m_SB_bam2, trans = plogis, select=2, ylim=c(0, 0.003), shift = coef(m_SB_bam2)[1], shade = TRUE, shade.col = "lightblue") #mat age
plot(m_SB_bam, trans = plogis, select=4, ylim=c(0.0005, 0.0015), shift = coef(m_SB_bam2)[1], shade = TRUE, shade.col = "lightblue") #month
##gam.check(m_SB_bam2)
concurvity(m_SB_bam2)## para s(GA_weeks) s(mat_age) s(month)
## worst 0.9818192 0.010299471 0.04741036 0.0005560260
## observed 0.9818192 0.005179840 0.01102309 0.0003709731
## estimate 0.9818192 0.006635969 0.03035472 0.0002292482 #adding birthweight
m_SB_bam3 = update(m_SB_bam2, . ~ . + s(BW))
summary(m_SB_bam3)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(month) + Urban3 + Language +
## mother_nationality_cat2 + mean_ssep2 + s(BW)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.780549 0.142433 -47.605 <2e-16 ***
## Urban3 0.029018 0.039126 0.742 0.4583
## LanguageFrench 0.022621 0.044259 0.511 0.6093
## LanguageItalian -0.065520 0.102628 -0.638 0.5232
## mother_nationality_cat2Africa 0.217818 0.093867 2.320 0.0203 *
## mother_nationality_cat2Asia -0.042533 0.102073 -0.417 0.6769
## mother_nationality_cat2Europe 0.066854 0.042722 1.565 0.1176
## mother_nationality_cat2Northern America -0.301652 0.288144 -1.047 0.2952
## mother_nationality_cat2Oceania -0.077984 0.724539 -0.108 0.9143
## mother_nationality_cat2Southern and Central America -0.045985 0.136210 -0.338 0.7357
## mean_ssep2 -0.005025 0.002344 -2.143 0.0321 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.293 7.998 382.174 <2e-16 ***
## s(mat_age) 2.791 3.530 8.635 0.0552 .
## s(month) 2.698 3.355 7.965 0.0527 .
## s(BW) 6.934 7.927 876.522 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.139 Deviance explained = 34%
## fREML = 1.5384e+06 Scale est. = 1 n = 1102617 plot(m_SB_bam3, trans = plogis, select=1, ylim=c(0, 0.04) ,shift = coef(m_SB_bam3)[1], shade = TRUE, shade.col = "lightblue") #GA
plot(m_SB_bam3, trans = plogis, select=2, ylim=c(0, 0.003), shift = coef(m_SB_bam3)[1], shade = TRUE, shade.col = "lightblue") #mat age
plot(m_SB_bam3, trans = plogis, select=3, ylim=c(0.0005, 0.0015), shift = coef(m_SB_bam3)[1], shade = TRUE, shade.col = "lightblue") #month
plot(m_SB_bam3, trans = plogis, select=4, ylim=c(0, 1), shift = coef(m_SB_bam3)[1], shade = TRUE, shade.col = "lightblue") #BW with whole CI
plot(m_SB_bam3, trans = plogis, select=4, ylim=c(0, 0.2), shift = coef(m_SB_bam3)[1], shade = TRUE, shade.col = "lightblue") #BW
gam.check(m_SB_bam3) 
##
## Method: fREML Optimizer: perf newton
## full convergence after 4 iterations.
## Gradient range [-1.054505e-06,0.00284269]
## (score 1538366 & scale 1).
## Hessian positive definite, eigenvalue range [0.3200282,2.64928].
## Model rank = 47 / 47
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 7.29 0.97 0.86
## s(mat_age) 9.00 2.79 0.97 0.93
## s(month) 9.00 2.70 0.97 0.32
## s(BW) 9.00 6.93 1.01 1.00 concurvity(m_SB_bam3)## para s(GA_weeks) s(mat_age) s(month) s(BW)
## worst 0.9818236 0.8787410 0.04886795 0.0005680043 0.8787533
## observed 0.9818236 0.4613031 0.01742516 0.0004029686 0.5160357
## estimate 0.9818236 0.4136715 0.03144524 0.0002441208 0.3329934 #to calculate the OR of stillbirth
or_gam(data=bevn_eco_in6, m_SB_bam2, pred= "mean_ssep2", percentage=20, slice = TRUE)## predictor value1 value2 perc1 perc2 oddsratio CI_low (2.5%) CI_high (97.5%)
## 1 mean_ssep2 23.63 36.24 0 20 0.93 0.89 0.97
## 2 mean_ssep2 36.24 48.85 20 40 0.93 0.90 0.96
## 3 mean_ssep2 48.85 61.46 40 60 0.93 0.92 0.93
## 4 mean_ssep2 61.46 74.07 60 80 0.93 0.94 0.91
## 5 mean_ssep2 74.07 86.68 80 100 0.93 0.96 0.89 AIC(m_SB_bam2, m_SB_bam3)## df AIC
## m_SB_bam2 25.85570 30465.16
## m_SB_bam3 32.52614 29626.622008-2010
m_SB_bam_2009<-bam(stillbirth~s(GA_weeks)+ s(mat_age) + s(birth_Y_M_num) + s(BW) + mean_ssep2 + Urban3 + Language + mother_nationality_cat2,family="binomial", data=bevn_eco_in6_2008_10)
summary(m_SB_bam_2009)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(BW) +
## mean_ssep2 + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.763035 0.320684 -21.089 <2e-16 ***
## mean_ssep2 -0.005448 0.005282 -1.031 0.3024
## Urban3 0.160198 0.087471 1.831 0.0670 .
## LanguageFrench -0.021483 0.101040 -0.213 0.8316
## LanguageItalian -0.009108 0.209587 -0.043 0.9653
## mother_nationality_cat2Africa -0.309746 0.263353 -1.176 0.2395
## mother_nationality_cat2Asia -0.638783 0.295570 -2.161 0.0307 *
## mother_nationality_cat2Europe 0.112526 0.094428 1.192 0.2334
## mother_nationality_cat2Northern America -0.786455 0.750256 -1.048 0.2945
## mother_nationality_cat2Oceania -8.028036 51.939069 -0.155 0.8772
## mother_nationality_cat2Southern and Central America -0.209763 0.317197 -0.661 0.5084
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 5.747 6.804 68.243 <2e-16 ***
## s(mat_age) 1.850 2.343 2.142 0.410
## s(birth_Y_M_num) 1.012 1.023 1.735 0.189
## s(BW) 5.315 6.356 260.124 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.161 Deviance explained = 35.2%
## fREML = 3.0468e+05 Scale est. = 1 n = 220464 #birth_Y_M_num LINEAR
m_SB_bam_2009_2 <- update(m_SB_bam_2009, . ~ . + birth_Y_M_num - s(birth_Y_M_num))
summary(m_SB_bam_2009_2)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(BW) + mean_ssep2 +
## Urban3 + Language + mother_nationality_cat2 + birth_Y_M_num
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.594385 0.343528 -19.196 <2e-16 ***
## mean_ssep2 -0.005447 0.005282 -1.031 0.3024
## Urban3 0.160201 0.087471 1.831 0.0670 .
## LanguageFrench -0.021452 0.101039 -0.212 0.8319
## LanguageItalian -0.009084 0.209586 -0.043 0.9654
## mother_nationality_cat2Africa -0.309716 0.263351 -1.176 0.2396
## mother_nationality_cat2Asia -0.638781 0.295569 -2.161 0.0307 *
## mother_nationality_cat2Europe 0.112513 0.094428 1.192 0.2334
## mother_nationality_cat2Northern America -0.786541 0.750262 -1.048 0.2945
## mother_nationality_cat2Oceania -8.028067 51.938880 -0.155 0.8772
## mother_nationality_cat2Southern and Central America -0.209726 0.317196 -0.661 0.5085
## birth_Y_M_num -0.005426 0.004091 -1.326 0.1847
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 5.747 6.804 68.243 <2e-16 ***
## s(mat_age) 1.850 2.343 2.142 0.41
## s(BW) 5.315 6.356 260.126 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.161 Deviance explained = 35.2%
## fREML = 3.0469e+05 Scale est. = 1 n = 220464 plot(m_SB_bam_2009_2, trans = plogis, select=1, ylim=c(0, 0.04) ,shift = coef(m_SB_bam_2009_2)[1], shade = TRUE, shade.col = "lightblue") #GA
plot(m_SB_bam_2009_2, trans = plogis, select=2, ylim=c(0, 0.004), shift = coef(m_SB_bam_2009_2)[1], shade = TRUE, shade.col = "lightblue") #mat age
plot(m_SB_bam_2009_2, trans = plogis, select=3, ylim=c(0, 0.4), shift = coef(m_SB_bam_2009_2)[1], shade = TRUE, shade.col = "lightblue") #BW 
2017-2020
#birth date as smoothed variable
m_SB_bam_covid<-bam(stillbirth~s(GA_weeks)+ s(mat_age) + s(BW) + s(birth_Y_M_num) + mean_ssep2 + Urban3 + Language + mother_nationality_cat2,family="binomial", data=bevn_eco_in6_2017_20)
summary(m_SB_bam_covid)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(BW) + s(birth_Y_M_num) +
## mean_ssep2 + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -7.038400 0.261776 -26.887 < 2e-16 ***
## mean_ssep2 -0.001171 0.004289 -0.273 0.784827
## Urban3 -0.135189 0.071765 -1.884 0.059594 .
## LanguageFrench 0.118289 0.079199 1.494 0.135290
## LanguageItalian -0.291285 0.230558 -1.263 0.206449
## mother_nationality_cat2Africa 0.526574 0.149991 3.511 0.000447 ***
## mother_nationality_cat2Asia 0.209034 0.167461 1.248 0.211938
## mother_nationality_cat2Europe 0.056427 0.078726 0.717 0.473525
## mother_nationality_cat2Northern America 0.346601 0.441780 0.785 0.432715
## mother_nationality_cat2Oceania -8.018614 46.335232 -0.173 0.862607
## mother_nationality_cat2Southern and Central America 0.076200 0.243805 0.313 0.754627
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.089 7.925 156.226 <2e-16 ***
## s(mat_age) 1.004 1.007 1.395 0.240
## s(BW) 5.465 6.658 198.222 <2e-16 ***
## s(birth_Y_M_num) 1.002 1.004 0.768 0.382
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.129 Deviance explained = 34.6%
## fREML = 4.6254e+05 Scale est. = 1 n = 333721 plot(m_SB_bam_covid, trans = plogis, select=1, ylim=c(0, 0.1) ,shift = coef(m_SB_bam_covid)[1], shade = TRUE, shade.col = "lightblue") #GA
plot(m_SB_bam_covid, trans = plogis, select=2, ylim=c(0, 0.003), shift = coef(m_SB_bam_covid)[1], shade = TRUE, shade.col = "lightblue") #mat age
plot(m_SB_bam_covid, trans = plogis, select=3, ylim=c(0, 0.4), shift = coef(m_SB_bam_covid)[1], shade = TRUE, shade.col = "lightblue") #BW 
plot(m_SB_bam_covid, trans = plogis, select=4, ylim=c(0.0005, 0.0015), shift = coef(m_SB_bam_covid)[1], shade = TRUE, shade.col = "lightblue") #birth_Y_M_num 
#birth date and mat age as linear variables
m_SB_bam_covid2<-bam(stillbirth~s(GA_weeks) + s(BW) + birth_Y_M_num + mat_age + mean_ssep2 + Urban3 + Language + mother_nationality_cat2,family="binomial", data=bevn_eco_in6_2017_20)
summary(m_SB_bam_covid2)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(BW) + birth_Y_M_num + mat_age +
## mean_ssep2 + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.977928 0.472786 -14.759 < 2e-16 ***
## birth_Y_M_num -0.002198 0.002508 -0.876 0.380943
## mat_age 0.008080 0.006852 1.179 0.238262
## mean_ssep2 -0.001171 0.004289 -0.273 0.784808
## Urban3 -0.135187 0.071764 -1.884 0.059597 .
## LanguageFrench 0.118292 0.079199 1.494 0.135279
## LanguageItalian -0.291274 0.230558 -1.263 0.206465
## mother_nationality_cat2Africa 0.526593 0.149990 3.511 0.000447 ***
## mother_nationality_cat2Asia 0.209045 0.167460 1.248 0.211913
## mother_nationality_cat2Europe 0.056438 0.078725 0.717 0.473438
## mother_nationality_cat2Northern America 0.346598 0.441780 0.785 0.432719
## mother_nationality_cat2Oceania -8.018585 46.335324 -0.173 0.862608
## mother_nationality_cat2Southern and Central America 0.076216 0.243804 0.313 0.754574
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.089 7.925 156.2 <2e-16 ***
## s(BW) 5.465 6.658 198.2 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.129 Deviance explained = 34.6%
## fREML = 4.6254e+05 Scale est. = 1 n = 333721 plot(m_SB_bam_covid2, trans = plogis, select=1, ylim=c(0, 0.1) ,shift = coef(m_SB_bam_covid2)[1], shade = TRUE, shade.col = "lightblue") #GA
plot(m_SB_bam_covid2, trans = plogis, select=2, ylim=c(0, 0.4), shift = coef(m_SB_bam_covid2)[1], shade = TRUE, shade.col = "lightblue") #BW 
Gestational age: continuous
#GA in days
m_GA_bam_days <- bam(GA_days ~ s(BW) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat2 + Language, data=bevn_eco_in7)
summary(m_GA_bam_days)##
## Family: gaussian
## Link function: identity
##
## Formula:
## GA_days ~ s(BW) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat2 +
## Language
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 276.04111 0.01752 15751.707 < 2e-16 ***
## parity_cat(1,2] -1.99202 0.01754 -113.596 < 2e-16 ***
## parity_cat(2,3] -2.49376 0.02687 -92.814 < 2e-16 ***
## parity_cat(3,20] -2.56052 0.04653 -55.033 < 2e-16 ***
## sex2 1.81258 0.01581 114.644 < 2e-16 ***
## Urban3 0.09180 0.01657 5.539 3.04e-08 ***
## mother_nationality_cat2Africa 0.60971 0.04764 12.798 < 2e-16 ***
## mother_nationality_cat2Asia -0.79664 0.04316 -18.460 < 2e-16 ***
## mother_nationality_cat2Europe -0.35922 0.01783 -20.147 < 2e-16 ***
## mother_nationality_cat2Northern America -0.28620 0.10175 -2.813 0.00491 **
## mother_nationality_cat2Oceania -0.32052 0.23325 -1.374 0.16939
## mother_nationality_cat2Southern and Central America -1.59307 0.05865 -27.161 < 2e-16 ***
## LanguageFrench 0.10488 0.01928 5.441 5.29e-08 ***
## LanguageItalian -0.20128 0.04297 -4.684 2.81e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(BW) 8.977 9.000 127518.85 <2e-16 ***
## s(mat_age) 8.044 8.641 347.60 <2e-16 ***
## s(birth_Y_M_num) 8.689 8.972 115.30 <2e-16 ***
## s(month) 6.958 8.043 12.76 <2e-16 ***
## s(mean_ssep2) 7.867 8.668 59.48 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.515 Deviance explained = 51.5%
## fREML = 3.8714e+06 Scale est. = 67.031 n = 1099328plot(m_GA_bam_days, all.terms = TRUE,pages=1, shade = TRUE, shade.col = "lightblue")
plot(m_GA_bam_days, select=1, shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_days)[1]) #BW
plot(m_GA_bam_days, select=2, ylim=c(265, 280),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_days)[1]) #mat age
plot(m_GA_bam_days, select=3, ylim=c(274, 278),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_days)[1]) #birth_Y_M_num
plot(m_GA_bam_days, select=4, ylim=c(275, 277),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_days)[1]) #month
plot(m_GA_bam_days, select=5, ylim=c(272, 280),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_days)[1]) #ssep
#checking that residuals are not determined
gam.check(m_GA_bam_days)
##
## Method: fREML Optimizer: perf newton
## full convergence after 13 iterations.
## Gradient range [-0.001200072,0.001194492]
## (score 3871436 & scale 67.03148).
## Hessian positive definite, eigenvalue range [1.310202,549654.5].
## Model rank = 59 / 59
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(BW) 9.00 8.98 1.02 0.87
## s(mat_age) 9.00 8.04 0.98 0.14
## s(birth_Y_M_num) 9.00 8.69 1.00 0.50
## s(month) 9.00 6.96 0.98 0.14
## s(mean_ssep2) 9.00 7.87 0.99 0.11 concurvity(m_GA_bam_days, full=TRUE)## para s(BW) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.8015268 0.04840718 0.12110433 0.03003088 0.028058758 0.17917903
## observed 0.8015268 0.01911969 0.07090416 0.01041702 0.001316444 0.04934839
## estimate 0.8015268 0.03133441 0.08072505 0.01100951 0.021033963 0.08954157#GA in weeks
m_GA_bam_weeks <- bam(GA_weeks ~ s(BW) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat2 + Language, data=bevn_eco_in7)
summary(m_GA_bam_weeks)##
## Family: gaussian
## Link function: identity
##
## Formula:
## GA_weeks ~ s(BW) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat2 +
## Language
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 39.434444 0.002504 15751.707 < 2e-16 ***
## parity_cat(1,2] -0.284574 0.002505 -113.596 < 2e-16 ***
## parity_cat(2,3] -0.356252 0.003838 -92.814 < 2e-16 ***
## parity_cat(3,20] -0.365788 0.006647 -55.033 < 2e-16 ***
## sex2 0.258940 0.002259 114.644 < 2e-16 ***
## Urban3 0.013115 0.002368 5.539 3.04e-08 ***
## mother_nationality_cat2Africa 0.087102 0.006806 12.798 < 2e-16 ***
## mother_nationality_cat2Asia -0.113806 0.006165 -18.460 < 2e-16 ***
## mother_nationality_cat2Europe -0.051318 0.002547 -20.147 < 2e-16 ***
## mother_nationality_cat2Northern America -0.040885 0.014536 -2.813 0.00491 **
## mother_nationality_cat2Oceania -0.045789 0.033322 -1.374 0.16939
## mother_nationality_cat2Southern and Central America -0.227581 0.008379 -27.161 < 2e-16 ***
## LanguageFrench 0.014983 0.002754 5.441 5.29e-08 ***
## LanguageItalian -0.028755 0.006139 -4.684 2.81e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(BW) 8.977 9.000 127518.85 <2e-16 ***
## s(mat_age) 8.044 8.641 347.60 <2e-16 ***
## s(birth_Y_M_num) 8.689 8.972 115.30 <2e-16 ***
## s(month) 6.958 8.043 12.76 <2e-16 ***
## s(mean_ssep2) 7.867 8.668 59.48 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.515 Deviance explained = 51.5%
## fREML = 1.7323e+06 Scale est. = 1.368 n = 1099328 plot(m_GA_bam_weeks, all.terms = TRUE,pages=1, shade = TRUE, shade.col = "lightblue")
plot(m_GA_bam_weeks, select=1,shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_weeks)[1])
plot(m_GA_bam_weeks, select=2, ylim=c(38, 40),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_weeks)[1])
plot(m_GA_bam_weeks, select=3, ylim=c(39.3, 39.5),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_weeks)[1])
plot(m_GA_bam_weeks, select=4, ylim=c(39.3, 39.5),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_weeks)[1])
plot(m_GA_bam_weeks, select=5, ylim=c(39.1, 39.9),shade = TRUE, shade.col = "lightblue", shift = coef(m_GA_bam_weeks)[1])
gam.check(m_GA_bam_weeks)
##
## Method: fREML Optimizer: perf newton
## full convergence after 13 iterations.
## Gradient range [-0.001200072,0.001194492]
## (score 1732279 & scale 1.367989).
## Hessian positive definite, eigenvalue range [1.310202,549654.5].
## Model rank = 59 / 59
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(BW) 9.00 8.98 1.02 0.86
## s(mat_age) 9.00 8.04 0.99 0.23
## s(birth_Y_M_num) 9.00 8.69 0.99 0.38
## s(month) 9.00 6.96 1.01 0.75
## s(mean_ssep2) 9.00 7.87 0.98 0.16 concurvity(m_GA_bam_weeks, full=TRUE)## para s(BW) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.8015268 0.04840718 0.12110433 0.03003088 0.028058758 0.17917903
## observed 0.8015268 0.01911969 0.07090416 0.01041702 0.001316444 0.04934839
## estimate 0.8015268 0.03133441 0.08072505 0.01100951 0.021033963 0.08954157#try with maternal age with k=12
m_GA_bam_1 = update(m_GA_bam_weeks, . ~ . - s(mat_age) + s(mat_age, k=13))
summary(m_GA_bam_1)##
## Family: gaussian
## Link function: identity
##
## Formula:
## GA_weeks ~ s(BW) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) +
## parity_cat + sex + Urban3 + mother_nationality_cat2 + Language +
## s(mat_age, k = 13)
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 39.434435 0.002504 15751.561 < 2e-16 ***
## parity_cat(1,2] -0.284571 0.002505 -113.595 < 2e-16 ***
## parity_cat(2,3] -0.356254 0.003838 -92.815 < 2e-16 ***
## parity_cat(3,20] -0.365793 0.006647 -55.034 < 2e-16 ***
## sex2 0.258942 0.002259 114.645 < 2e-16 ***
## Urban3 0.013116 0.002368 5.540 3.03e-08 ***
## mother_nationality_cat2Africa 0.087110 0.006806 12.799 < 2e-16 ***
## mother_nationality_cat2Asia -0.113808 0.006165 -18.460 < 2e-16 ***
## mother_nationality_cat2Europe -0.051316 0.002547 -20.147 < 2e-16 ***
## mother_nationality_cat2Northern America -0.040873 0.014536 -2.812 0.00493 **
## mother_nationality_cat2Oceania -0.045810 0.033322 -1.375 0.16920
## mother_nationality_cat2Southern and Central America -0.227581 0.008379 -27.161 < 2e-16 ***
## LanguageFrench 0.014980 0.002754 5.440 5.33e-08 ***
## LanguageItalian -0.028756 0.006139 -4.684 2.81e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(BW) 8.977 9.000 127518.99 <2e-16 ***
## s(birth_Y_M_num) 8.688 8.972 115.30 <2e-16 ***
## s(month) 6.958 8.043 12.76 <2e-16 ***
## s(mean_ssep2) 7.868 8.668 59.48 <2e-16 ***
## s(mat_age) 9.472 10.434 287.93 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.515 Deviance explained = 51.5%
## fREML = 1.7323e+06 Scale est. = 1.368 n = 1099328gam.check(m_GA_bam_1)
##
## Method: fREML Optimizer: perf newton
## full convergence after 13 iterations.
## Gradient range [-0.001574849,0.001567499]
## (score 1732280 & scale 1.367989).
## Hessian positive definite, eigenvalue range [1.310621,549654.5].
## Model rank = 62 / 62
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(BW) 9.00 8.98 1.01 0.76
## s(birth_Y_M_num) 9.00 8.69 0.99 0.37
## s(month) 9.00 6.96 1.01 0.69
## s(mean_ssep2) 9.00 7.87 1.01 0.83
## s(mat_age) 12.00 9.47 1.03 0.99#pvalue not significant, k-index for month and birthdate is < 1.
#try with mat age k=12 and month k=10
m_BW_bam4_2 <- bam(BW ~ s(GA_weeks) + s(mat_age, k=13) + s(birth_Y_M_num) + s(month, k=11) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language + mother_nationality_cat2, data=bevn_eco_in7)
summary(m_BW_bam4_2)##
## Family: gaussian
## Link function: identity
##
## Formula:
## BW ~ s(GA_weeks) + s(mat_age, k = 13) + s(birth_Y_M_num) + s(month,
## k = 11) + s(mean_ssep2) + parity_cat + sex + Urban3 + Language +
## mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3335.2684 0.8361 3988.889 < 2e-16 ***
## parity_cat(1,2] 127.4692 0.8345 152.758 < 2e-16 ***
## parity_cat(2,3] 173.6232 1.2772 135.937 < 2e-16 ***
## parity_cat(3,20] 199.6911 2.2149 90.157 < 2e-16 ***
## sex2 -147.0104 0.7458 -197.129 < 2e-16 ***
## Urban3 0.6595 0.7906 0.834 0.404
## LanguageFrench -48.0926 0.9182 -52.374 < 2e-16 ***
## LanguageItalian -63.1901 2.0493 -30.835 < 2e-16 ***
## mother_nationality_cat2Africa 11.0235 2.2735 4.849 1.24e-06 ***
## mother_nationality_cat2Asia -17.9212 2.0596 -8.702 < 2e-16 ***
## mother_nationality_cat2Europe 32.5734 0.8501 38.316 < 2e-16 ***
## mother_nationality_cat2Northern America 65.0139 4.8546 13.392 < 2e-16 ***
## mother_nationality_cat2Oceania 63.4734 11.1288 5.704 1.17e-08 ***
## mother_nationality_cat2Southern and Central America 47.0408 2.7991 16.806 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(GA_weeks) 8.828 8.988 73215.198 < 2e-16 ***
## s(mat_age) 6.568 7.836 14.641 < 2e-16 ***
## s(birth_Y_M_num) 8.285 8.859 39.006 < 2e-16 ***
## s(month) 5.978 7.248 3.072 0.00277 **
## s(mean_ssep2) 7.634 8.534 43.676 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.401 Deviance explained = 40.1%
## fREML = 8.1205e+06 Scale est. = 1.5259e+05 n = 1099328 gam.check(m_BW_bam4_2)
##
## Method: fREML Optimizer: perf newton
## full convergence after 11 iterations.
## Gradient range [-0.0009357942,0.0007208003]
## (score 8120452 & scale 152595).
## Hessian positive definite, eigenvalue range [1.648423,549654.5].
## Model rank = 63 / 63
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(GA_weeks) 9.00 8.83 1.01 0.665
## s(mat_age) 12.00 6.57 0.98 0.065 .
## s(birth_Y_M_num) 9.00 8.28 1.00 0.365
## s(month) 10.00 5.98 0.99 0.270
## s(mean_ssep2) 9.00 7.63 1.02 0.905
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1concurvity(m_BW_bam4_2, full=TRUE) ## para s(GA_weeks) s(mat_age) s(birth_Y_M_num) s(month) s(mean_ssep2)
## worst 0.8015963 0.017211066 0.12307439 0.030165524 0.028070793 0.17908212
## observed 0.8015963 0.006386136 0.04255397 0.007228201 0.002163185 0.10239354
## estimate 0.8015963 0.009760255 0.07758667 0.011469079 0.019775195 0.08954325 #from the "worst" case, the values are not high (high would be >0.8, usually)
plot(m_BW_bam4_2, select=3, ylim=c(-15, 20), xlim=c(144, 168),shade = TRUE, shade.col = "lightblue")
#with mat nationality category 3
m_GA_bam_2 <- bam(GA_weeks ~ s(BW) + s(mat_age) + s(birth_Y_M_num) + s(month) + s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat3 + Language, data=bevn_eco_in7)
summary(m_GA_bam_2)##
## Family: gaussian
## Link function: identity
##
## Formula:
## GA_weeks ~ s(BW) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + parity_cat + sex + Urban3 + mother_nationality_cat3 +
## Language
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 39.451420 0.002509 15723.245 < 2e-16 ***
## parity_cat(1,2] -0.280711 0.002513 -111.719 < 2e-16 ***
## parity_cat(2,3] -0.345646 0.003856 -89.631 < 2e-16 ***
## parity_cat(3,20] -0.349566 0.006656 -52.518 < 2e-16 ***
## sex2 0.259523 0.002259 114.883 < 2e-16 ***
## Urban3 0.015323 0.002369 6.468 9.91e-11 ***
## mother_nationality_cat3Central Europe -0.059664 0.004532 -13.166 < 2e-16 ***
## mother_nationality_cat3Eastern Europe -0.049213 0.011336 -4.341 1.42e-05 ***
## mother_nationality_cat3Northern Europe -0.098784 0.017882 -5.524 3.31e-08 ***
## mother_nationality_cat3other -0.069801 0.003994 -17.478 < 2e-16 ***
## mother_nationality_cat3SE Europe -0.109024 0.004266 -25.559 < 2e-16 ***
## mother_nationality_cat3Southern Europe -0.062939 0.005616 -11.207 < 2e-16 ***
## mother_nationality_cat3SW Europe 0.036876 0.005412 6.813 9.55e-12 ***
## mother_nationality_cat3Western Europe 0.015800 0.006592 2.397 0.0165 *
## LanguageFrench 0.002541 0.002841 0.894 0.3711
## LanguageItalian -0.031802 0.006183 -5.144 2.69e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(BW) 8.977 9.000 127489.25 <2e-16 ***
## s(mat_age) 8.088 8.671 377.13 <2e-16 ***
## s(birth_Y_M_num) 8.687 8.971 120.50 <2e-16 ***
## s(month) 6.994 8.073 12.82 <2e-16 ***
## s(mean_ssep2) 7.839 8.653 62.58 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.514 Deviance explained = 51.4%
## fREML = 1.7325e+06 Scale est. = 1.3684 n = 1099368Gestational age: PTB
m_SB_bam<-bam(stillbirth~s(GA_weeks)+ s(mat_age) + s(birth_Y_M_num) + s(month) +s(mean_ssep2) + Urban3 + Language + mother_nationality_cat2,family="binomial", data=bevn_eco_in6)
summary(m_SB_bam)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(birth_Y_M_num) + s(month) +
## s(mean_ssep2) + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.96844 0.04200 -165.927 <2e-16 ***
## Urban3 0.03816 0.03892 0.980 0.3269
## LanguageFrench 0.08108 0.04404 1.841 0.0656 .
## LanguageItalian -0.02042 0.10239 -0.199 0.8420
## mother_nationality_cat2Africa 0.22796 0.09350 2.438 0.0148 *
## mother_nationality_cat2Asia -0.03662 0.10164 -0.360 0.7187
## mother_nationality_cat2Europe 0.03212 0.04252 0.755 0.4500
## mother_nationality_cat2Northern America -0.34845 0.28509 -1.222 0.2216
## mother_nationality_cat2Oceania -0.17323 0.72070 -0.240 0.8100
## mother_nationality_cat2Southern and Central America -0.09875 0.13512 -0.731 0.4649
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.589 8.143 16734.682 < 2e-16 ***
## s(mat_age) 3.227 4.052 13.034 0.01214 *
## s(birth_Y_M_num) 1.013 1.026 0.010 0.95230
## s(month) 2.627 3.267 7.214 0.07017 .
## s(mean_ssep2) 1.008 1.016 6.865 0.00925 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.133 Deviance explained = 32.1%
## fREML = 1.5655e+06 Scale est. = 1 n = 1102617plot(m_SB_bam)
plot(m_SB_bam, pages = 1, trans = plogis)
#adding intercept
plot(m_SB_bam, pages = 1, trans = plogis,
shift = coef(m_SB_bam)[1])
plot(m_SB_bam, trans = plogis, select=1,,shift = coef(m_SB_bam)[1], shade = TRUE, shade.col = "lightblue")
plot(m_SB_bam, trans = plogis, select=2, ylim=c(0, 0.003), shift = coef(m_SB_bam)[1], shade = TRUE, shade.col = "lightblue")
plot(m_SB_bam, trans = plogis, select=3, ylim=c(0.0005, 0.0015), shift = coef(m_SB_bam)[1], shade = TRUE, shade.col = "lightblue")
plot(m_SB_bam, trans = plogis, select=4, ylim=c(0.0005, 0.0015), shift = coef(m_SB_bam)[1], shade = TRUE, shade.col = "lightblue")
plot(m_SB_bam, trans = plogis, select=5, ylim=c(0.0005, 0.0015),shift = coef(m_SB_bam)[1], shade = TRUE, shade.col = "lightblue")
# second model_ without date of birth (Only month), and putting mean_ssep2 as linear term
m_SB_bam2 = update(m_SB_bam, . ~ . - s(mean_ssep2) + mean_ssep2 - s(birth_Y_M_num))
summary(m_SB_bam2)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks) + s(mat_age) + s(month) + Urban3 + Language +
## mother_nationality_cat2 + mean_ssep2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.614331 0.141044 -46.896 < 2e-16 ***
## Urban3 0.038204 0.038917 0.982 0.32626
## LanguageFrench 0.081151 0.044036 1.843 0.06535 .
## LanguageItalian -0.020692 0.102355 -0.202 0.83979
## mother_nationality_cat2Africa 0.228259 0.093470 2.442 0.01460 *
## mother_nationality_cat2Asia -0.036356 0.101616 -0.358 0.72051
## mother_nationality_cat2Europe 0.032281 0.042502 0.760 0.44754
## mother_nationality_cat2Northern America -0.348623 0.285086 -1.223 0.22138
## mother_nationality_cat2Oceania -0.173444 0.720702 -0.241 0.80982
## mother_nationality_cat2Southern and Central America -0.098668 0.135115 -0.730 0.46523
## mean_ssep2 -0.006102 0.002331 -2.617 0.00886 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 7.589 8.143 16736.825 <2e-16 ***
## s(mat_age) 3.224 4.049 13.026 0.0121 *
## s(month) 2.627 3.268 7.209 0.0703 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.133 Deviance explained = 32.1%
## fREML = 1.5655e+06 Scale est. = 1 n = 1102617 plot(m_SB_bam2, pages = 1, trans = plogis)
#adding intercept
plot(m_SB_bam2, pages = 1, trans = plogis,
shift = coef(m_SB_bam)[1])
plot(m_SB_bam2, trans = plogis, select=1,,shift = coef(m_SB_bam2)[1], shade = TRUE, shade.col = "lightblue")
plot(m_SB_bam2, trans = plogis, select=2, ylim=c(0, 0.003), shift = coef(m_SB_bam2)[1], shade = TRUE, shade.col = "lightblue")
plot(m_SB_bam, trans = plogis, select=4, ylim=c(0.0005, 0.0015), shift = coef(m_SB_bam2)[1], shade = TRUE, shade.col = "lightblue")
##gam.check(m_SB_bam2)
concurvity(m_SB_bam2)## para s(GA_weeks) s(mat_age) s(month)
## worst 0.9818192 0.010299471 0.04741036 0.0005560260
## observed 0.9818192 0.005179840 0.01102309 0.0003709731
## estimate 0.9818192 0.006635969 0.03035472 0.0002292482 sum_m_SB_bam2<- summary(m_SB_bam2)
m_SB_bam2$p.table## NULLp_table <- data.frame(sum_m_SB_bam2$p.table)
p_table <- within(p_table, {lci <- Estimate - qnorm(0.975) * Std..Error
uci <- Estimate + qnorm(0.975) * Std..Error})
p_table## Estimate Std..Error z.value Pr...z..
## (Intercept) -6.614331055 0.141043601 -46.8956481 0.000000000
## Urban3 0.038203756 0.038917302 0.9816651 0.326264879
## LanguageFrench 0.081151266 0.044036189 1.8428312 0.065353655
## LanguageItalian -0.020692424 0.102354783 -0.2021637 0.839788729
## mother_nationality_cat2Africa 0.228259387 0.093469753 2.4420669 0.014603441
## mother_nationality_cat2Asia -0.036355524 0.101616309 -0.3577725 0.720513552
## mother_nationality_cat2Europe 0.032280907 0.042501608 0.7595220 0.447540339
## mother_nationality_cat2Northern America -0.348622684 0.285085797 -1.2228694 0.221379045
## mother_nationality_cat2Oceania -0.173444374 0.720701638 -0.2406604 0.809818302
## mother_nationality_cat2Southern and Central America -0.098668374 0.135114921 -0.7302552 0.465234204
## mean_ssep2 -0.006101804 0.002331376 -2.6172546 0.008864019
## uci lci
## (Intercept) -6.337890677 -6.89077143
## Urban3 0.114480267 -0.03807275
## LanguageFrench 0.167460611 -0.00515808
## LanguageItalian 0.179919263 -0.22130411
## mother_nationality_cat2Africa 0.411456736 0.04506204
## mother_nationality_cat2Asia 0.162808782 -0.23551983
## mother_nationality_cat2Europe 0.115582527 -0.05102071
## mother_nationality_cat2Northern America 0.210135210 -0.90738058
## mother_nationality_cat2Oceania 1.239104881 -1.58599363
## mother_nationality_cat2Southern and Central America 0.166152004 -0.36348875
## mean_ssep2 -0.001532391 -0.01067122 m_SB_bam3<-bam(stillbirth~s(GA_weeks, k=15)+ s(mat_age, k=15) + s(month) +mean_ssep2 + Urban3 + Language + mother_nationality_cat2,family="binomial", data=bevn_eco_in6)
summary(m_SB_bam3)##
## Family: binomial
## Link function: logit
##
## Formula:
## stillbirth ~ s(GA_weeks, k = 15) + s(mat_age, k = 15) + s(month) +
## mean_ssep2 + Urban3 + Language + mother_nationality_cat2
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.616927 0.141079 -46.902 <2e-16 ***
## mean_ssep2 -0.006057 0.002332 -2.597 0.0094 **
## Urban3 0.039152 0.038931 1.006 0.3146
## LanguageFrench 0.079560 0.044053 1.806 0.0709 .
## LanguageItalian -0.019312 0.102361 -0.189 0.8504
## mother_nationality_cat2Africa 0.229611 0.093465 2.457 0.0140 *
## mother_nationality_cat2Asia -0.036058 0.101609 -0.355 0.7227
## mother_nationality_cat2Europe 0.032767 0.042512 0.771 0.4408
## mother_nationality_cat2Northern America -0.353583 0.285739 -1.237 0.2159
## mother_nationality_cat2Oceania -0.167973 0.719951 -0.233 0.8155
## mother_nationality_cat2Southern and Central America -0.100768 0.135293 -0.745 0.4564
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(GA_weeks) 11.011 12.410 16739.723 <2e-16 ***
## s(mat_age) 3.263 4.156 13.117 0.0131 *
## s(month) 2.645 3.289 7.325 0.0675 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.133 Deviance explained = 32.2%
## fREML = 1.5619e+06 Scale est. = 1 n = 1102617 ##gam.check(m_SB_bam3)
AIC(m_SB_bam2, m_SB_bam3)## df AIC
## m_SB_bam2 25.85570 30465.16
## m_SB_bam3 29.84624 30448.53